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Thermodynamics of clusterized matter is studied in the framework of statistical models with 
non-interacting cluster degrees of freedom. At variance with the analytical Fisher model, exact 
Metropolis simulation results indicate that the transition from homogeneous to clusterized matter 
lies along the p = po axis at all temperatures and the limiting point of the phase diagram is not 
a critical point even if the surface energy vanishes at this point. Sensitivity of the inferred phase 
diagram to the employed statistical framework in the case of finite systems is discussed by considering 
the grand-canonical and constant-pressure canonical ensembles. A Wigner-Seitz formalism in which 
the fragment charge is neutralized by an uniform electron distribution allows to build the phase 
diagram of neutron star matter as predicted by the cluster model. 
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I. INTRODUCTION 



Droplet models constitute a suitable tool to study structural ordering transitions. Freezing, condensation and 
glassiness together with metallic alloys phase transitions illustrate both the phenomena diversity and the multitude 
of domains which benefit from them |l| . 

At sub-atomic scale, clusterized nuclear matter is expected to be present in the cores of exploding supernovae 
@, 0, 0, 0] and in the inner crust of neutron stars 0, 0, 0, 0, E3, 11 , 12, HI where the interplay between surface and 
Coulomb energy is responsible for the occurrence of exotic pasta phases at densities close to the saturation density of 
nuclear matter. 

It is conjectured that clusterized nuclear matter is formed also at the break-up stage of multifragmentation reactions, 
i.e. the simultaneous multi-particle decay of excited matter produced in energetic heavy ion reactions. It is this 
scenario that persuaded Finn et at 14 L 15[ to address nuclear multifragmentation is terms of a critical evaporation 
according to the Fisher cluster theory 16| as early as the '80s. The approach was daring the more as a parallelism 
among a neutral infinite system and an electrically-charged small one is far from being trivial. Nonetheless, more and 
more refined multifragmentation data continued to be scrutinized for criticality signals over more than two decades 
13, e, m, m [2l|, | §2]. Most of the results tend to confirm the criticality scenario and point-out that a first-order phase 
transition occurs as well 23|, 24, 25, 26|, 27],[28|]. However a real proof of criticality when analyzing a finite system like 
a fragmenting nucleus would need finite size scaling techniques, and these techniques are not experimentally accessible 
as the size of the system cannot be varied freely in nuclear fragmentation. 

From a theoretical point of view, several statistical models with cluster degrees of freedom [2^, |3(J [HI, H3] have 
been proposed and proved able to describe a wealth of experimental data in great detail. This opens the possibility 
of addressin g th e issue of criticality from a theoretical point of view. However past studies have only concerned finite 
systems 3^, 34, 35 1 , while the thermodynamic limit is needed to observe non- ambiguous signals of criticality. In this 
respect these results are not conclusive. 

The aim of the present paper is to investigate the phase diagram and the critical behavior of an infinite neutral 
clusterized system. To reach this goal we undertake the cluster model but, at variance with its standard analytical 
representation commonly known as the Fisher model we solve it exactly via a Metropolis simulation. This 
approach will allow us to see to what extent the analytical approximation stands and, secondly, to investigate how 
the different physical ingredients act on the thermodynamical properties. More precisely, we shall focus on effects of 
free volume correction, fragment compressibility, internal excitation, translational degrees of freedom and finite size 
effects. We shall additionally show that the cluster model can present first and second order phase transitions but the 
limiting point of the phase diagram is never a critical point even if the surface energy is made vanish at this point. 
This implies that either the cluster model in the different versions presented in the literature is not appropriate to 
describe nuclear fragmentation, or that this latter phenomenon cannot be interpreted as critical. 

The paper is organized as follows. Sec. II discusses the thermal and phase properties of nuclear matter as obtained 
by an exact Metropolis Monte-Carlo model with respect to the predictions of the analytical Fisher model. Sec. Ill 
analyzes the sensitivity of small systems to the considered thermodynamical framework, while the phase diagram of 
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neutron star matter in the Wigner-Seitz approximation is presented in Sec. IV. Conclusions are drawn in Sec. V. 

II. THE CLUSTER MODEL 

It is theoretically well settled that uniform nuclear matter is instable to density fluctuations of finite wavelength 
in a wide region of densities and temperatures, leading to spontaneous cluster formation 36|, 3?||. The same is true 



for stellar matter [381 . 139t I40H. and the physical properties at the crust-core transition point arc of special interest 
for astrophysical purposes [41|. A description of the thermodynamical properties of these systems can be obtained 
introducing a statistical model with cluster degrees of freedom. Such a n ap proach was proposed by Fisher in the late 
sixties to describe the process of condensation close to a critical point [161 ]. 

Under the equilibrium assumption, the statistical approach of a system broken into TV pieces reduces the physical 
problem to the estimation of the number of microscopic states compatible with the thermodynamical macroscopic 
constraints. For instance, in the simple case of a 1-component fluid, if one ignores the internal excitation, position, 
momenta, spins and parities of the clusters, a configuration will be exclusively defined by fragment sizes, 

k:{A u ...,A Nk }. (1) 

Including isotopic, translational and internal degrees of freedom, as it is usually done to treat the particular case 
of nuclear multifragmentation, leads to: 

k : {A\, Zi,ei,fi,pi, . . . , A Nk , Z Nk , e Nk , r Nk , p Nk } , (2) 

where Ai, Zi, e,, r*, and pi stand for the mass, charge, internal energy, center-of-mass position and momentum of each 
fragment. 

Assuming that it is possible to write down a tractable expression for the statistical weight of the configuration Wk 
in the statistical ensemble suitable for the phenomenon under study, one can express the characteristic partition sum 
Z and any ensemble-averaged observable < X > as a function of the statistical weight of the configuration Wk , 



Z = J2 W ^ ( 3 ) 

(ft) 



and, respectively, 



<x> = v w ■ (4) 

Z^(fc) W k 

Many improvements and sophistications can and have to be added to describe the nuclear multifragmentation prob- 
lem [32| or the stellar equation of state [1, [f| , but in this section we shall only focus on one-component clusterized 
systems in the grandcanonical ensemble. This will allow for a direct comparison between the exact solution of the clus- 
ter model within Monte-Carlo methods and an approximation allowing analytical expressions at the thermodynamic 
limit. 

A. The analytical Fisher model 

In the grand-canonical ensemble, the generic expression of the partition function of an infinite one-component 
system composed of independent clusters reads, 

= ]T cxp [-/? (e£> - ^ot)] , (5) 

(fc) 

with, 

00 

E^ = J2^e(A); (6) 

A=l 

00 

42 = E^ (?) 

A=l 
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where is the multiplicity of fragments of size A in the (fc)-class of events, e(A) is the cluster energy functional, 
and for the moment we ignore the degeneracy of the different energy states for each fragment size. The interest of the 
Fisher cluster model [3] is that the partition sum can be expressed as a product of individual cluster partition sums, 

oo 

zp,^ = n z wMi (8) 

A=l 

where 

zp t pM) = X] ex P[-/^(e(^) 

nA 

= 1 ± exp [-(3(e(A) - pA)\ , (9) 

and the + (— ) accounts for Fermi (Bose) statistics of the clusters states. If we consider the classical Boltzmann limit 
of a huge number of available states, ha <C 1 we get, 

up,Pfi(A) = In zp t pn(A) « cxp [-f3(e(A) - fiA)] , (10) 

leading to, 



oo oo n A / a\ 

= n e <"> 



A=l n A =0 

or, equivalcntly, 



hxZpjn = "p^A), ( 12 ) 



A=l 



with ujp^p tl (A) given by Eq. |[H 

Taking now into account that the clusters may be excited and have translational degrees of freedom e n (A 1 p) = 
e n (A) +p 2 /2Am, Eq. ([12]) becomes, 

In Z ,p, « £ f [ f^Ej29n(A)e X p[-(3(e n (A,^-fxA)] (13) 



A=l 



(rp \ 3/2 oo 
^A^exph^^-M)], (W) 



where g n (j4) accounts for the degeneracy of internal state (n), the internal entropy is given by exps(j4, e) = 
J2 n 9n{A)S(e — e n (A)), and we have made a saddle point approximation, 

deexp [s(A, e) - 0e] « fp(A), (15) 

where the cluster free energy at temperature T = is given by fa = (e(A))p — T(s(A))p. 

It is important to stress that the original Fisher model (Ref. [l(|) employs the saddle point approximation Eq. 
(TT5")) . but it does not include translational and internal degrees of freedom. For consistency these latter should be 
accounted directly in the partition sum Eq. (JS]), see Eq. (J3J). If this is done the key expression Eq. © is not valid 
any more, and this is why we include these extra degrees of freedom only at the level of the single-cluster partition 
sum ujp t p ll {A) 1 leading to the approximation Eq. (I13|) . 

For nuclear (fermionic) clusters both the average cluster energy and entropy can be evaluated in the low temperature 
Fermi gas approximation 

(e(A))p = -B(A) + a AT 2 ; (16) 
(s(A))p = (2a r + a ^- 1 / 3 /(T)) A -r\n A, (17) 

and to have an analytically solvable model, the fragment binding energy B = —e\ can be simply parameterized as: 

B(A) = a v A — a s A 2 / 3 . (18) 
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The surface term a s A 2 / 3 f(T) in Eq. (JXTj) effectively accounts for the entropy increase at finite temperature due to 
surface excitations. In particular if we choose for f(T) a smooth function such that f(T) = (3 for f3 < (3 C , this naturally 
produces a vanishing surface free energy at the critical point (3c = Tq 1 . Concerning the Fisher topologic factor r In A, 
it insures a fragment fractal dimension at the critical point. The physical meaning of this term, as pointed out in Ref. 

is related to the treatment of the center of mass degree of freedom close to a critical point. In such a situation, 
the available states for the center of mass motion have to be reduced respect to those of a free particle, and this can 
be accounted by reducing the exponent associated to the translational energy factor cx A 3 / 2 in Eq. (fT4|) . We have 
taken tq = r — 3/2 = 2.2 in order to recover the fractal dimension of the LG universality class. A more sophisticated 



temperature dependent expression for r has been proposed in Ref. 42 1. 

The specific form of Eqs. (fT7|) and (fT6|) insures that the model is critical at T = Tc, H = fic = ~ a i> — o-oTq. Indeed, 
computing the average multiplicity of clusters of size A, 

i \ dlnZpfj, _ 

{nA) ^ = -am = = 

= v(^) 3/2 A 3 / 2 e X p[-(3(MA)-^A)}. (19) 



One can see that the partition sum of the Fisher model is simply given by the sum of all the possible cluster 
multiplicities: 

oo 

In Zpfin = ( n ^)/3.^- ( 20 ) 

A=l 

At T = Tc, fJ. = He the fragment multiplicity obeys a power law, 



(mT c \ 



3/2 



(nA) 0c , 0ofi o = V\^J A~^. (21) 

Since the different thermodynamic observables can be expressed as successive derivatives of the partition sum, this 
means that all thermodynamic quantities will exhibit a power law behavior at the approach of T = Tc, /i = \ic- bi 
particular the critical exponents 7, a may be calculated numerically out of the total particle number, 



X cx \T c -T\-\ 



T T da 2 



(22) 



and, respectively, total energy fluctuation, 



C v oc \T C -T\ 

_ 1 d 2 In Zp !/3ll 



T 2 T 2 d/3 2 



(23) 



It is easy to see that any value for the topologic exponent 2 < tq < 3 will produce a diverging susceptibility x an d 
a finite critical density, 



3/2 

r - 2 V 2ttH 2 ) 



The critical pressure is also immediately calculated as, 

Pc~T c -Tc[^) ^— T (25) 

This formulation together with its canonical and microcanonical counterparts have been extensively used to stud 



nuclear multifragmentation and most of the interesting features of its associated thermodynamics [43j, |44 



to s tuay 



phase transition, critical behavior, finite size, isospin asymmetry and Coulomb effects, often with contradictory results. 
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B. The Metropolis Monte-Carlo approach 



A complementary way to solve our problem is by employing a Metropolis Monte-Carlo investigation of the con- 
figuration space. While having the obvious disadvantage of being analytically non-tractable, this method has the 
advantage of being solvable without approximations for any finite number of particles. Let us start again from the 
grancanonical partition sum of the cluster model for a finite system of volume V: 



Zp,fo = 2^n7\J J^N k exp r /? [ E tot - V A tot )\ , (26) 



with 



E%> =Y,e{A i: p i ,n i ); (27) 



Aol = Y. A » ( 28 ) 

i=l 

where again the cluster energy functional consists of kinetic energy and internal levels: e{A 1 p 1 n) = e n (A) + p 2 /2Am. 
Calculating the phase space integral Eq. (|26[) simplifies to 

(k) 



(k) k i=l V 7 



expl-PifdiAJ-vLAi)], (29) 



where we can take the same free energy functional Eqs. (fTol [T7|) as for the Fisher model in order to compare the 
two approaches. In principle inter-fragment interactions should be taken into account in the form of an excluded 
volume [33|. This effect is not accounted in the original Fisher model because it is a size-dependent effect which 
becomes negligible in the thermodynamic limit of interest here, as we will discuss in the next section. Therefore we 
will consider point-like particles also in the Monte-Carlo version of the model, Eq. 
In the thermodynamic limit we can write 



N k 

i-En=nE- (30) 

(fc) i=l A=ln A =0 



This means that, in the absence of the factorials 1/iVfc! in Eq. (|2"9"|) and l/n^! in Eq. the analytical and 

Monte-Carlo model would be equivalent at the thermodynamic limit. Because of the factorials, this equivalence is 
violated: considering the exact counting of all the possible configurations including their translational and internal 
degrees of freedom breaks the independence of the different cluster sizes Eq. © even if we keep the Fisher hypothesis 
of non-interacting point-like clusters. Eq. (|20[) has then to be considered as an approximation of the exact partition 
sum of the cluster model in the thermodynamic limit. The thermodynamic properties of the cluster model should 
then be investigated solving numerically Eq. (f2"9")l for a finite size system big enough that the independence of the 
system size is achieved, or that converging quantities can be extracted through finite size scaling techniques. 

The spanning of the available configurations in Eq. (|29p is done by proposing successively trial configurations which 
are accepted according to the detailed balance principle: 

N k \P(k- > k')W k Ak = N k AP(k'- > k)W k ,Ak', (31) 

see Ref. [32[ for details. Direct access to full thermodynamical information of each statistical configuration allows one 



to trace phase coexistence from the bimodal structure of the different probability density distributions [471 . l48j and, 
moreover, to identify to which phase each configuration belongs. 

In general and by construction, the double peak structure of the probability distribution originates from the convex 
intruder in the surface of the associated density of states as a function of the order parameter of the underlying first- 



order phase transition 47] ■ We expect particle and energy density as possible order parameters. In the grand-canonical 
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case under consideration, the probability distribution in the presence of a first order phase transition should then be 
double-peaked with respect to the extensive obscrvables A to t and E to t conjugated to the Lagrange multipliers (3 and 
a = /3[a. 

To allow for a straightforward comparison with the case of nuclear matter, in this Section, the following set of 
parameters will be used: (a v =16 MeV, a s = 16 MeV, a =l/16 MeV -1 , T c =16 MeV). To test the results stability to 
model parameters, two different expressions have been employed for f(T): 



(i) f(T) = T- 

(ii) f(T) = T- 



1 

Mi 



((T£-T 2 )/(T£4 
(1 - T/T c ) 6(T C 



T 2 )) 1 ' 25 Q(Tc 
-T)]. 



T) 



and 



C. The phase diagram 

Let us consider a system confined into the finite volume V and denote by Aq the maximum number of particles 
corresponding to normal density, pq = Aq/V. 

Fig. [1] illustrates the probability distributions of the total number of particles (left panel) and total energy (right 
panel) for a system characterized by 1^=14476.4 fm 3 (A =2000) at different values of (/3, a), together with the phase 
diagram (stars) in temperature - total density (left) and temperature - total energy (right) representations. To allow 
superposition of the above plots, the normalization factors of the probability distributions have been arbitrarily altered 
such that the peak heights correspond to the coexistence temperature represented on the left Y-axis. To estimate 
finite size effects, the phase diagram of a system two times smaller (1^=7238.2 fm 3 and Ao=1000) is represented with 
open squares. 

Given that Atot and Etot are correlated, at the thermodynamic limit the phase-coexistence information provided 
by the bimodal structure of the these distributions should be the same. Indeed, the phase-coexistence temperature 
- defined by the temperature value for which the peaks belonging to the dense and diluted phases have the same 
height if a is kept constant, differs by less than 0.1 MeV when one considers Y{A to t) or Y(E to t) distributions. This 
is coherent with a one-dimensional order parameter transition like ordinary liquid-gas. 

It is interesting to observe that correct thermodynamic information is obtained also from the probability distri- 
butions of the largest cluster A max in each event (middle panel) as in fragmentation transitions. This is due to the 
strong correlation of A max with the total energy deposit (49l. IBol. l5l| . Even more interesting, plotted as a function of 
(1 — A max /Ao), the phase diagram built out of A rnax bimodality sits perfectly on the top of the one built analyzing 
the bimodality of Y{A to t) distributions. Thus, we confirm that A max is a good order parameter in the cluster model. 
This is an encouraging result for the experimental search of the finite nuclei phase diagram out of the measured A max 
distributions (27I. Eq. 

As one may see in Fig. [TJ the shape of the phase diagram does not depend on the system size. The gas-like branch 
approaches the critical point of the analytical Fisher model discussed in the previous section. However, the limiting 
point characterized by Tu m =16 MeV and au rn =-2, where by construction the fragment size distribution is an exact 
power low, lies on the line p = po, which is very different from the value p c ~ 0.07488 fm -3 obtained in the analytical 
Fisher model Eq. |24|) . The same result has been obtained in Ref. [4f|, where the phase diagram of a 1-component 
clusterized system was calculated analytically out of the Gibbs equations. This result is not physical and is due to 
the fact that the cluster model considers all nuclear clusters fully incompressible, including the asymptotic limit of 
nuclear matter. 

Within this modclisation, the shape of the phase diagram and in particular its limiting temperature are independent 
of the system size for any Aq > 100. This means that the thermodynamic limit may be addressed directly from finite 
size systems calculations without performing finite size scaling. Finite size scaling implies a power law for the size 
dependence of the limiting temperature if this latter corresponds to a critical point, according to 

(T Hm (L) - T c ) L- 1 ^ , (32) 

where L is the system linear dimension, Tc = lim/,_ >00 Tu m (L), and v is the critical exponent related to the divergence 
of the correlation length at the critical point, £ cx |1 — T/Tc\~" '. The observed independence of Tu m (L) from the 
system size indicates that the correlation length of the infinite system does not diverge at Tc, which means that the 
limiting point is not a critical point. This again is at variance with the results of the analytical Fisher model, and in 
contradiction with most expectations on nuclear multifragmentation based on statistical clusters models. 

Before leaving this section, we mention that the shape of the phase diagram proves insensitive to the functional 
dependence according to which the surface energy term vanishes while approaching Tc- The solid stars and open 
squares in Fig. [T] correspond to the case in which f(T) has the form (i). 
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FIG. 1: (color online) Probability distributions of the total number of particles (left), largest fragment in each event (middle) 
and total energy (right) for a one-component cluster model with standard Fisher fragment definition in the grand-canonical 
ensemble characterized by V= 14476.4 fm 3 (,4o=2000) and (T, (ji/T)) = (8.3, -2.40), (10.35, -2.20) and (13.19, -2.05), where T is 
expressed in MeV (Metropolis simulation results). The normalization factors are chosen such that the peak heights correspond 
to the coexistence temperature. The solid star symbols mark the phase diagram in temperature - density (left and middle 
panel) and temperature - total energy representations as obtained out of the bimodality of corresponding distributions. Open 
squares mark the phase diagram of a system half the size, (V=7238.2 fm 3 and A o =1000). 



D. Criticality 

The question whether or not this second-order transition point is indeed a critical point may be also answered cal- 
culating different critical exponents and checking whether they do verify, together with the input r, the corresponding 
equality relations [l6[ , 

7 - (3 - T) = i-^r. (33) 



(2 - a) (r - 1) 2 - a 

and. 

The full thermodynamic characterization of the system allows the calculation of (3 and 5 out of the curvature of the 
coexistence curve, 

and, respectively, the shape of the critical isotherm, 

{pc - P) \t c <* (Pc ~ p) 5 \t c , (36) 
which is equivalent through (pc — p) oc pc(pc — p) to, 

(Mc - P) \t c oc {p c - p) S \t c ■ (37) 

The absence of finite size effects in the phase diagram allows direct extraction of (3 out of the curvature of the 
coexistence curve, Eq. (pI5j) . The results of the Metropolis simulation for the one-component clusterized system in 
the grandcanonical ensemble are plotted in panel (a) of Fig. [2] with solid stars along with the predictions of the 
analytical model, represented with open circles. Over the explored (p/T) domain the coexistence curve p CO ex exhibits 
the expected power law behavior as a function of T in the two models, with slightly different values for the exponents. 
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By contrast with this, when approaching the critical point the critical isotherm deviates from the expected linear 
behavior and the magnitude of these deviations increases with the system size, as illustrated in the panel (b) of Fig. 
[5] where full and open stars correspond to the predictions of the grandcanonical one-component Metropolis model 
with ^0=2000 and 4000. This behavior can be understood taking into account that the calculation of p demands the 
evaluation of the average number of particles, a quantity which, like any statistical average, is naturally affected by 
the finite size of the system. The slope of the linear region 5=5.2 is, nevertheless, fairly close to the predictions of the 
analytical model such that one may conclude that the hyperscaling equation (f34f holds. 

7 and a may be directly extracted out of the fluctuations of the system total energy and total number of particles 
via Eqs. (j22|) and (f23|) along the critical a=-2 path. 

The results of the Fisher analytical model shown in panels (c,d) fulfill again perfectly the hyperscaling relations 
Eq. (|33p . as expected from the discussion of the previous section. Since fluctuations should diverge at the approach 
of a critical point, the Monte-Carlo results obtained for a finite system must obviously saturate at a finite value. 
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FIG. 2: (color online) Extraction of critical exponents corresponding to a one-component clusterized system within a grand- 
canonical ensemble out of the thermodynamical properties. Results corresponding to the analytical model are plotted with 
open circles, while stars stand for the results of the Metropolis simulation. In this last case, the maximum number of particles 
allowed by the fixed volume is Ao = 2000 (solid stars) and, respectively, Ao — 4000 (open stars). Panel a: /3 extracted out of 
the coexistence curve p — T; panel b: <5 extracted out of the critical isotherm; panel c: 7 extracted out of total particle number 
fluctuations; panel d: a extracted out of total energy fluctuations. Cases (c) and (d) consider the critical {p/T) path. 
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TABLE I: Critical exponents for a one component clusterized system within a grandcanonical ensemble as obtained from a 
Metropolis simulation in comparison with results of the analytical model and the expectations from the Liquid-Gas universality 
class. 



Critical exponent 


Calculation method 


Monte Carlo Analytical 


LG 


a 


a%/T 2 oc \T-T c \~ a 


0.34 


0.2 


0.1 


7 


a\/T cc IT-Tcl" 7 


2.8 


1.2 


1.2 


P 


(pc-pcoe*) oc (T c -T) 


0.25 


0.33 


0.33 


8 


(pc -p)\t c « (pc - p) 5 \t c 


5.2 


5 


5 


7/r 


mi oc (T c — T)~ 7 


2.8 




1.2 


Pfr 


A max oc (T c -T) fl 


0.33 




0.33 


Whhr 


A or ^ /3/7 

Jt-max ^c *->2 


0.16 


0.275 


0.275 



This can be seen in panels c and d, where, irrespective of the system size, one may distinguish two regimes for the 
Monte-Carlo results: in the temperature domain 0.1 < 1 — T/Tc < 0.7, <j\/T and a 2 E /T 2 are linear with a non-null 
slope, while for T approaching Tc, the two curves get flat. The comparative analyze of A o =2000 and 4000 results 
shows that neither the widths of the two intervals, neither the value of the non-null slope depend on the system size. 
The first observation is a first indication that the system might not be critical. Indeed, if the system were critical, 
by increasing the system size, the linearity domain close to the critical point would increase. This is obviously not 
the case here. The second observation is in line with the lack of system size dependence in the phase diagram and 
allows one to extract the critical exponents without performing finite size scaling. The values a = 0.34 and 7 = 2.8 
obtained over the region where ln(l — p coe x/pc) vs. hi(l — T/Tc) is linear are so different from the expected values 
a = 0.2,7 = 1.2 extracted from Eq. (j3"3")l that it is difficult to believe that this can be a finite size effect. 

The ensemble of these results show that criticality is violated in the cluster model even if the surface tension vanishes 
at the limiting point of the first-orderphase transition. 

In percolation theory over infinite and small lattices [5^ , critical exponents may be extracted from fc-moments 
of the cluster distributions, 

mi = ^2 sn ( s ) (P ~ Pc) , (38) 



m 2 



E" 5 



s n(s) oc\p- pc\~ 



(39) 



where p, the probability that a lattice site is occupied, is the control parameter. When dealing with finite systems at 
T < Tc, the largest cluster assimilated with the liquid phase must be excluded from summations [52| . 

In nuclear multifragmentation, and Eqs. (|38|) and (|39p have been employed in the equivalent forms 17, lj| 



(Tc-T) 13 , for T <T C 



(40) 



m 2 cx (T c - T)-\ for T < T c , (41) 

and, 

A ma x CX ^2^1 (42) 

where S2 = rn 2 / 'mi, and often the distance to the critical point T — Tc has been estimated from the fragment 
multiplicity M. Eq. (|42p is particularly easy to handle as it does not require data sorting according to T , a quantity 
very difficult to infer experimentally. 

Fig. [3] investigates to what extent the parallelism between percolation on lattice and clusterized matter holds. 
The evolution of the second moment of the cluster distribution and the largest fragment in each event as a function 
of (1 — T/Tc) for T < Tc are plotted along the critical ev=-2 path (panels a and b). After eliminating the region 
where divergences are rounded-off by finite size effects, one gets "ff r = 2.8 and /3/ r =0.33. The result for 7 is in 
excellent agreement with the fluctuation estimation, and again very different from the value 7 = 1.2 one should get 
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to have thermodynamically consistent exponents. The slight disagreement between the values may be surprising 
considering that we have observed that A max is an excellent order parameter in the cluster model. However it is 
interesting to note that, at variance with all other procedures here employed, Eq. (|4T)|) proves sensitive to the way 
in which fragment surface energy vanishes while approaching the limiting temperature. For instance for f(T) = 
T- 1 [1 - (1 - T/T c ) Q{T C - T)], Eq. (gDJ) leads to fr =O.25. The predictions of Eq. (02J are illustrated in panel c 
for 12 McV < T < 24 MeV along a=-2. Liquid-like and gas-like branches of the Campi plot give (0/j)f r =O.16. This 
number is not compatible with criticality which would require (0H) = 0.275, but it is in striking agreement with the 
values reported in the experimental multifragmentation literature 13, 2^|- This result is surprising since deviations 
might have been expected in data due to the Coulomb interaction. 

Before leaving this section let us mention that, if we would only dispose of cluster observables as in the experimental 
case of multifragmentation, the results for {f3/-f)f r , 7/r an d 0f r from Eqs. (|4*0l . ([ST]) and (|4"2"jl . would make us 
erroneously conclude that the exponents are thermodynamically consistent and the system is critical. The computation 
of thermodynamical quantities as in Figure [2] above is necessary to realize the violation of criticality of the model. 
Therefore there is no contradiction between the fact that the cluster model is not critical and the fact that this model 
is so successful in reproducing experimental data which have been qualified as critical. 

The values of the above discussed critical exponents are summarized in Table I together with the employed equations. 



E. Temperature dependent bulk energy 

In this study we have always considered completely incompressible point-like fragments. In this section we focus 
on what happens to the thermodynamics of a one-component clusterized system if, in addition to a s , also a v vanishes 
while approaching Tq- From the physical point of view, no microscopic calculation suggests that this can be the case, 
however this simplified calculation may give a glimpse on what could happen if compressibility effects were realistically 
introduced, leading to a reduced bulk energy at finite temperature. The phase diagram obtained as described in Sec. 
IIC is represented in Fig. [¥J with solid circles. The same symbols are used in Figs. [5] and [5] to depict the evolution 
of different thermodynamical and fragment observables out of which critical exponents are calculated as described in 
Sec. IID. 

The first result is that, when plotted in T — p, the phase diagram superimposes perfectly on the top of the one 
corresponding to standard fragment definition. This is however not valid in T — E representation, because of the 
different fragment energetics and, moreover, (p/T)c=0. The stability of the T— p phase coexistence borders suggests 
that the critical exponents and 6 extracted out of the curvature of phase coexistence and shape of the critical 
isotherm have the same values as before. The confirmation is given in Fig. [5j panels a and b. Panel d shows that 
neither the divergence of ki(a%/T 2 ) vs. ln(l — T/Tc) is affected. By contrast, ln(er^/T) vs. ln(l — T/Tc) gets a linear 
dependence over a broader region and 7=1.9. Hyperscaling relations are still violated and critical behavior does not 
hold. Fig. [6] shows that critical exponents extracted out of fragment properties are in reasonable agreement with the 
values calculated from the thermodynamic behavior. In this case, f(T) = T^ 1 [1 — (1 — T/Tc) @(Tc — T)]. 



F. Free volume correction 



Working under the point-like cluster hypothesis, as we did so far, the system's thermodynamics is entirely dictated 
by fragment energetics. More sophisticated cluster-based statistical models include, in the spirit of Van der Waals 
theory of fluids, the geometrical effect of clusters finite size which effectively act as a repulsive fragment-fragment 
interaction for incompressible clusters. 

To understand which is the role the excluded volume effect plays in the present model, let us consider what happens 
if each cluster forbids for all the other clusters a volume equal with its own volume. Mathematically, this means that 
the factor V Nk entering in the Eq. (|29[) should be replaced with 

V Nk -> V(V - Vx)(V -V x - V 2 )...(V -Vi-V 2 - ... - Vzvfc-i), (43) 

where Vi is the intrinsec volume of fragment i, here considered as being spherical, Vi = Air/Sv^Ai. 

Fig. [7] shows the phase diagram of a Aq = 2000 1-component system obtained as described in Sec. IIC. We note 
that, for the first time, the liquid branch of the phase coexistence region is situated at p < po and, because of geometric 
effects, it bends at T > 8 MeV. This second behavior obviously prevents criticality studies. 
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III. APPLICATION TO FINITE SYSTEMS - SENSITIVITY TO THE STATISTICAL ENSEMBLE 



As we have already mentioned, Fisher cluster models have been intensively used to address multifragmentation 
[43L 44, |45| and the thermodynamic behavior of finite systems has been constantly compared to the one of their 
infinite counterparts. The long-range non-saturating Coulomb interaction, which has to be included when describing 
finite fragmenting nuclei, prevents addressing the thermodynamic limit and will not be considered here. Interference 
with isospin-asymmctry effects is also avoided by exclusively considering isospin-symmetric systems. Even if we ignore 
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FIG. 3: (color online) Extraction of critical exponents 7 (Eq. (|39)l . top panel), /3 (Eq. (|40)l . middle panel) and j3/*y (Eq. (|42)l . 
bottom panel) out of fragment properties along the critical (/i/T)=-2 path in the case of a one-component clusterized system 
within a grand-canonical ensemble. The maximum number of particle allowed by the fixed volume is Ao — 2000. 
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FIG. 4: The phase diagram in temperature - density (left and middle panel) and temperature - total energy representations as 
obtained out of the bimodality of corresponding distributions for a system characterized by V=14476.4 fm 3 (v4o=2000) within 
the grand-canonical ensemble. At variance with standard fragment definition, here we assume a temperature dependence of 
the bulk energy. 



these effects, realistic applications to fragmenting finite nuclear systems differ from the model we have presented in 
Section II also because they take into account the excluded volume correction to the translational motion, do not 
generally consider the topologic Fisher factor and treat the internal cluster level density realistically without employing 
the low-temperature approximation eq. (|17p . 

These differences suggest that the multifragmcntation phase diagram may be different from the one plotted in Fig. 

m 

To see whether it is the case, the phase diagrams of different charge-neutral isospin-symmetric systems {Aq— 200, 
500, 3600) obtained as described in Section IIC are depicted in Fig. [5] with open circles, open stars and, respectively, 
solid circles. Standard nuclear multifragmentation fragment definition is employed [32j. Parameters of fragment 
binding energy have values close to the one used in Section II, a„=15.4941 MeV, a s =17.9439 MeV, and we have 



additionally taken for the symmetry energy ai = 1.7826 [54(. Nuclear level density is taken as 



= 12a i/4 £ 5/4 cx p( 2 ^) exp(-e/r), (44) 

with a = 0.1144 + 0.098A 2 / 3 MeV" 1 and r=9 MeV 

As one may notice, the liquid border of the coexistence zone is stable at p = pa, but the temperature of the 
second-order transition point is situated well below Tp=16 MeV and, moreover, presents a system size dependence. 
The shift is essentially due to the factor A 3 / 2 in the partition sum coming from the translational motion which does 
not vanish even at Tq- In this case, it is clear that not only critical behavior is lost, but also the power law shape of 
fragment size distributions in the vicinity of the phase diagram limiting point. This result is in qualitative agreement 
with the conclusions of Ref. [43| , where the translational energy factor in the fragment partition sum is responsible 
for the diminish of the second-order transition temperature from 16 MeV to 7 MeV and the lost of linear shape of 
fragment size distributions in log-log scale. 

The system size dependence of Tn m is due to the inclusion of inter-fragment interactions by means of an excluded 
volume in the spirit of the Van der Waals gas and persists up to Aq w 4000. 

Indeed, if we consider that only the space not occupied by the clusters is available for the center-of-mass motion, the 
factor in Eq. has to be replaced by Eq. (gSJ). Assuming for simplicity that Vi = V 2 = ■ ■ ■ = V Nk = (V/n)/N k , 
which is reasonable only for relatively symmetric fragment partitions or high multiplicities, Eq. (|43[) may be written 
as V Nk £ with, 

iVfc-l / . «. 

{ -n('-^). <«> 
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FIG. 5: (color online) The same as in Fig. [2] for the case in which bulk energy depends on temperature. Results of Metropolis 
calculation (solid circles) correspond to Ao=2000. 



where n = V/Vq, see Ref. |32j for details. Since, 



lim £ 

V— >oo,Vb— *oo 



1. 



(46) 



we can see that the excluded volume effect becomes negligible only for very large systems. This confirms that the 
point-particle approximation is correct only when dealing with the thermodynamic limit. 

Even more important, the thermodynamics of a finite system depends on the statistical framework in which it 
is investigated. A constant pressure canonical ensemble may be used to trace p — T phase diagrams out of the 
bimodality of the volume and to immediately compare them with the phase diagrams obtained grandcanonically out 
of the bimodality of the total number of particles. The results corresponding to two systems (^4=200 and 1000) are 
depicted in Fig. [9] with open and solid circles in different representations. As above, the system size dependence 
of Tn m is the consequence of excluded volume corrections. Contrary to the grandcanonical results, in this case the 
liquid border of the phase coexistence region extends over the whole density domain and the coexistence region at 
low temperatures shrinks to the vicinity of p = 0. The limiting point is characterized by p = po. As the size of the 
system increases, the constant-pressure canonical phase diagram approaches the grand-canonical one to which it is 
identical in the thermodynamic limit. The convergence is nevertheless very slow. 
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FIG. 6: (color online) The same as in Fig. [3] for the case in which bulk energy depends on temperature. Results of Metropolis 
calculation are represented with solid circles and Ao=2000. 



IV. APPLICATION TO INFINITE SYSTEMS: NEUTRON STARS 



Clusterized nuclear matter in the temperature and density domain of a liquid-gas phase transition is pred icted to 
be produced also in the core of supernovae and in the inner crust of neutron stars @, 0, S, 0, E3, 11, 12, EH- In 
this case, the system's net charge neutrality is insured by the ultra-relativistic gas of electrons in which the nuclear 
fragments are embedded. 

In the free neutron regime (4.3 x 10 11 g/cm 3 < p < 2.5 x 10 14 g/cm 3 ) the matter is presently figured-out as a 
lattice of nuclei immersed in a nucleon and electron gas and theoretically treated as a succession of non-interacting 
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FIG. 7: (color online) The phase diagram in temperature - density (left panel) and temperature - total energy representations 
as obtained out of the bimodality of corresponding distributions for a system characterized by V— 14476.4 mi 3 (^0=2000) 
within the grand-canonical ensemble. A free volume correction is included according to Eq. (|43l) . 



Wigner-Seitz cells 51 [II [57 1. 

If wc adopt this scenario, the thermodynamics of the whole star crust may be investigated analyzing only one such 
a cell and represents a perfect application of the cluster model. The effect of the Coulomb energy is that of modifying 
the energy of cluster configuration, and of suppressing the critical behavior of nuclear matter Because of that 
last point, we have not included the Fisher topologic factor — Tin A in the cluster entropy. Concerning the modified 
energetics, this effect is taken into account by altering the statistical weight of a configuration Eq. (f2"9"| according to, 



exp(-pV c (k))V N <°l 



mAiT 
2irh 2 



3/2 



(47) 



where the weight associated to a given cluster (Ai, Zi) is 

TXnwp^.^A^ZA = (a v Ai - a s (1 - Tf(T)) A 2 / 3 ) (l - a,(A t - 2Z i ) 2 /A i ) 
- e(Ai, Zi) + n n (Ai - Zi) + fipZi, 



(48) 



and nuclear fragments are described as in Sec. III. The Coulomb energy may be calculated for each configuration as 

El, 



v c {k) = Y,Up)—/ 

i=i b r oA 



1/3' 



with 



c(p) 



1 


-) - 




(Pe_\ 




KPopJ 




KPop) 



(49) 



(50) 



accounting for the screening effect of electrons. po p = Z/Apo denotes the proton density inside the nuclei and p e is 
the electron density. 

Fig. [TOl illustrates the phase diagram of a Wigner-Seitz cell of volume F=28952.91 fm 3 (Ao=4000) as obtained from 
the bimodality of total particle number distributions in a grandcanonical ensemble. 

One may directly compare the obtained phase diagram with the one corresponding to a system of similar size 
where Coulomb is switched-off (solid symbols in Fig. [8|). We can see that the shape of the phase diagram does not 
change. This would not necessarily be true any more if the weight of the translational motion oc (VAT) 3 / 2Nk at high 
temperature was reduced through a temperature dependent reduction factor as in the Lattimcr-Swesty equation of 
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FIG. 8: The phase diagram of the cluster model including excluded volume effects according to Eq. (|45[l as obtained out of 
the bimodal behavior of the total number of particles within the grandcanonical ensemble. Left side: a(— /3fi) as a function of 
density. Right side: potential energy per particle, a and density as a function of temperature. Empty circles, empty stars and 
solid circles correspond to charge-neutral isospin-symmetric systems with volume V=1447.64 fm 3 (Aq = 200), V=3617. 28 fm 3 
(A o =500) and, respectively, V = 26057.62 fm 3 (A = 3600). 



state [7[. However the prescription for such a term is completely phcnomenological, and we will not employ it here. 
A clear effect is seen on the limiting temperature, which is diminished by several MeV, in agreement with results 
obtained in Ref. [34[ for small systems. An opposite trend has been obtained within a classical microscopic model for 
Coulomb frustration and the understanding of this discrepancy requires further analyses. 

As observed in section II, the fact that the high density transition line lies along p = po is not physical and is due to 
the incompressibility assumption of the dense phase which is here represented as a single uncharged nucleus of infinite 
size. This line cannot therefore be interpreted as a prediction for the crust-core transition density as a function of 
temperature. Moreover, at densities close to saturation, exotic pasta phases are known to appear 0, Hoi , 
and the cluster energy functional should be modified to account for these non-spherical geometries. However in the 
low density region this model should give a good approximation of the thermodynamics of proto-neutron stars and 
supernova matter @, 0, Si • The evaluation of the equation of state with the inclusion of gammas and neutrinos is in 
progress. 



17 




0.2 

CO 

£0.15 
> 

o 0.1 



q_0.05 




O (200,100) 
• (1000,500) 



' Ma 



la; 



25 50 75 100 

v/v 



0.2 
0.15 
0.1 
0.05 

1 

0.75 









Q • 




# 




p * 




a / 













:: • (c) 




i I i i i I i i i 




6 8 10 

T (MeV) 
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V. CONCLUSIONS 



In conclusion, the thermodynamics of clusterized matter was investigated in connection with a liquid-gas phase 
transition in nuclear matter. The phase diagram presents first and second order phase transitions and may be derived 
out of the grand-canonical bimodality structure of total particle number, total energy and largest cluster in each 
event. Results of an exact Metropolis simulation are systematically compared with predictions of the analytical 
Fisher model which by construction is critical. At variance with the analytical counterpart, the predictions of the 
exact model indicate that the liquid phase is characterized by p = p$ and that vanishing fragment surface energy is 
not a sufficient condition to obtain criticality. 

The clusterized system model is applied to charge-neutral isospin-symmetric finite systems and special attention is 
given to the way in which the system thermodynamics depends on the employed statistical framework. Clusterized 
systems with net charge neutrality may be used also to mimic star matter. A first schematic calculation of the 
phase diagram in the Wigner-Seitz cell approximation is built and it is shown that within this approximation charge 
fluctuations lead to the decrease of the limiting temperature. 
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